Slow rupture of frictional interfaces 
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The failure of frictional interfaces and the spatio-temporal structures that accompany it are central 
to a wide range of geophysical, physical and engineering systems. Recent geophysical and laboratory 
observations indicated that interfacial failure can be mediated by slow slip rupture phenomena which 
are distinct from ordinary, earthquake-like, fast rupture. These discoveries have influenced the way 
we think about frictional motion, yet the nature and properties of slow rupture are not completely 
understood. We show that slow rupture is an intrinsic and robust property of simple non-monotonic 
rate-and-state friction laws. It is associated with a new velocity scale c m i n , intrinsically determined 
by the friction law, below which steady state rupture cannot propagate. We further show that 
rupture can occur in a continuum of states, spanning a wide range of velocities from c m in to elastic 
wave-speeds, and predict different properties for slow rupture and common fast rupture. Our results 
are in qualitative agreement with recent high-resolution laboratory measurements and may provide 
a theoretical framework for understanding slow rupture phenomena along frictional interfaces. 



I. INTRODUCTION 

Understanding the dynamic processes that govern in- 
terfacial failure and frictional sliding, e.g. a large earth- 
quake along a natural fault, remains a major scientific 
challenge. Quite recently, several geophysical and high- 
resolution real-time laboratory observations have pointed 
to the possibility that stress releasing interfacial slip can 
be mediated by the propagation of rupture fronts whose 
velocity is much smaller than elastic wave-speeds [IHEj]. 
These slow rupture fronts typically radiate seismic en- 
ergy at significantly reduced rates and hence are less de- 
structive than ordinary earthquakes. Furthermore, the 
propagation of slow rupture fronts was observed to pre- 
cede major slip events, possibly during the nucleation 
and early stages of large scale rupture, and hence raised 
some interest in their potential predictive roles. 

The nature and properties of these slow rupture fronts, 
and in particular their propagation velocity, are still 
not fully understood. The experiments of 0, Q clearly 
demonstrate the existence of a minimal propagation ve- 
locity below which no fronts are observed. To the best 
of our knowledge, no theoretical understanding of this 
minimal velocity is currently available. 

Frictional phenomena is commonly described using 
phenomenological rate-and-state friction models 
Two possible mechanisms for generating slow rupture 
events were invoked in this framework. The first involves 
a non-monotonic dependence of the steady state fric- 
tional resistance on slip velocity , while the second 
involves spatial variation of friction parameters and the 
accompanying heterogeneous stress distribution [l2l [l3j . 
The former mechanism is an intrinsic property of the 
friction law, while the latter mechanism is an extrin- 
sic one. The high-resolution laboratory measurements 
of [J,!, performed on a quasi-2D spatially homoge- 
neous system, suggests that the second mechanism is not 
necessary for the existence of slow rupture. 

In this Letter we show that slow rupture naturally 
emerges in the framework of spatially homogeneous rate- 



and-state friction models. Our analysis is based on 
a friction model that includes an elastic response at 
small shear stresses and a transition to slip above a 
threshold stress, following the "minimal" microscopic 
model of and its continuum generalization |T(J] . The 
model exhibits a crossover from velocity-weakening fric- 
tion at small slip rates to velocity-strengthening friction 
at higher slip rates, which we argue to be a generic fea- 
ture of friction. 

The existence of a minimal rupture front velocity c m i n , 
which is determined by the friction law and is indepen- 
dent of elastic wave-speeds, is predicted analytically in a 
quasi- ID limit. We show that there exists a continuum 
of rupture fronts with velocities ranging between c m ;„ 
and elastic wave-speeds, in qualitative agreement with 
recent laboratory measurements [l| and consistent with 
field observations [4| . We further show that slow rupture 
is significantly less spatially localized than ordinary fast 
rupture. These predictions are corroborated by explicit 
calculations for a glassy polymer that was extensively 
used in laboratory experiments [l], [H, [3 [I?} , demonstrat- 
ing the existence of slow rupture which is well-separated 
from ordinary fast rupture. We propose that these re- 
sults are possibly relevant for slow/silent earthquakes in 
geological contexts. 



II. A RATE-AND-STATE FRICTION MODEL 



Our goal here is to further develo p a nd incorporate 
the ideas discussed recently in [l5|, [H, [l8[ into a realistic 
rate-and-state friction model. We consider a spatially ex- 
tended interface between macroscopic bodies in frictional 
contact. The interface is composed of an ensemble of con- 
tact asperities whose total area A r is much smaller than 
the nominal contact area A n and which exert a shear 
stress t that resists sliding motion. We decompose r 
into an elastic part, emerging from the elastic deforma- 
tions of contact asperities which are characterized by a 
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coarse-grained stress T el , and a viscous part r ms 

t = r el + r ms = r el + r, v*A sgn(w) log ^1 + , (l) 

where 77 is a viscous-friction coefficient, v is the slip ve- 
locity (rate), v* is a reference velocity scale and A — 
A r /A n <C 1 is the normalized real contact area. The 
viscous-stress r ms , which increases with v and scales with 
A, is usually associated with activated rate processes at 
asperity contacts Q. The 1 inside the log ensures a reg- 
ular behavior in the limit v — > 0, but otherwise does not 
play a crucial role in what follows. 

The next step is writing down a dynamic evolution 
equation for r el . r el is stored at contact asperities at 
a rate determined by v and is proportional to both the 
intcrfacial elastic modulus fiQ and A. It is released as con- 
tact asperities are destroyed after slipping over a typical 
distance D (of the order of the size of a contact), when 
the asperity-level stress surpasses a yield-like threshold 
t c . This physical picture is mathematically captured by 
writing [16| 
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where h is the asperities' effective height. Note that the 
coarse-grained stress r is enhanced by a factor A^ 1 ^> 1 
at the asperities level and that the geometric nature of 
elastic stress relaxation, emerging from the multi-contact 
nature of the interface, is captured by the introduction of 
a spatial length D, that is interpreted as the typical size 
of a contact as in conventional rate-and-state models 
8]. The appearance of a step function non-linearity #(•) 
is an outcome of the basic notion of a static threshold 
for sliding motion. The evolution law in Eq. ([2]) features 
a reversible elastic response at small shear stresses, r ~ 
T el = noAu/h, where u is the shear displacement. This 
elastic response is usually not included in friction models 
(but see [19|, 12(1 ). even though it was directly measured 
experimentally [17| ■ We will later use these experimental 
data to constrain parameters such as fio/h. 

To proceed, we write the normalized contact area A in 
terms of a state variable d>, 



A{<t>,cx) = A {(j) [l + 6Iog(l+0/0o)] 



(3) 



Here a is the (compressive) normal stress and Aq(<j) = 
a /&h, where oh is the hardness. The evolution of the 
contact area is phenomenologically determined from two 
robust observations: logarithmic aging in the absence 
of irreversible slip and a logarithmic velocity-weakening 
steady-state behavior [8j . This is commonly captured by 
Dieterich's evolution law [a], extended here by stipulat- 
ing that the transition from the aging regime (v = 0) to 
the sliding regime (v ^ 0) is controlled by the same step 
function nonlinearity as in equation ([2]), yielding 



where 4> is interpreted as the "geometric age" of the con- 
tacts. Equations (HJ)- (J3J) determine the evolution of r(i), 
i.e. constitute our proposed friction law. 

Before we proceed we note a very important feature 
of rate-and-state friction models, which is not specific to 
the present model. In the absence of persistent sliding, 
v = 0, we have <f> = t and the contact area ages loga- 
rithmically Aoc 1 4- Mog(l-K/0o), as is widely observed 
Q. The latter form suggests that the logarithmic law 
is cutoff at short timcscales, smaller than O , as was 
directly confirmed experimentally [lil |2lj . This very 
same short timescales cutoff manifests itself also under 
persistent sliding, v ^ 0, for which we have <fr = D/\v\ 
and A oc 1 + 61og[l + D/((f>o\v\)]. In this case, A satu- 
rates at a finite value above a typical slip rate of order 
D/4>q and the fixed-point of Eq. ^ , r el oc A, becomes v- 
independent as well. As a consequence, r, which usually 
exhibits a velocity- weakening behavior at small v, be- 
comes velocity-strengthening as the viscous-friction term 
in Eq. (|TJ) takes over. Thus, rate-and-state friction mod- 
els generically predict a non-monotonic dependence of the 
steady state sliding friction on the slip rate, an observa- 
tion that has been largely overlooked in the literature 
(but see p| [13, HH, (13) an d that will play an important 
role below. 



III. STEADY STATE RUPTURE FRONTS 

At this point we are interested in steady state propa- 
gating front solutions corresponding to our friction law. 
Such propagating solutions exist in multi-stable systems 
in which one homogeneous (space independent) solution 
invades another one, giving rise to non-trivial spatiotem- 
poral structures. The spatially homogeneous solutions of 
equations ([T])-([4]), as a function of a driving stress r d , are 
shown in Fig. [T^,. A branch of elastic (static) solutions 
exists at v = 0, where aging effects are neglected, i.e. we 
assume that ipo = 61og(l + t/4>o) is roughly constant for 
the timescales relevant for front propagation. A branch 
of steady sliding solutions with v > takes the form 

/ = — ^ /o + "log (l + 4) +f 3h & f 1 + ~T~) . ( 5 ) 
where 

a = r)v*/a H , = fi Db/a H h , f = /3/b, (6) 

and / is the steady sliding friction coefficient. Note that 
we neglected a term of order log 2 in Eq. (JSJ) - As discussed 
above, steady sliding friction is indeed non-monotonic 
(when a</3)} friction is velocity- weakening for v*<S^v<^ 
D/cj)o and velocity-strengthening for v 3> D/<j>Q, with a 
minimum at 
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At the minimum, we define the friction stress as r m = 
T S s{v m ). Figs. [Tp-d present experimental data for three 
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FIG. 1: (a) A schematic sketch of the homogeneous solutions 
of t(v). The elastic (static) branch is shown as a vertical line 
at v = 0, which ends at the threshold t c A. The non-monotonic 
steady sliding curve is shown at v > 0. The slip velocity at the 
minimum is denoted by v m and the corresponding steady state 
friction stress by r m . The horizontal dashed lines denote dif- 
ferent levels of the driving stress r d with respect to r m and the 
corresponding homogeneous solutions (filled symbols), (b ) 
The steady sliding friction coefficient f(v) for PMMA [23 ] ■ 
The solid line is a fit to Eq. {5). This fit was highly con- 
strained by additional independent data sets, see Sect. IIVI 
(c) Af(v) (defined with respect to some reference level) for 
granite with <r = 5 MPa [H]. (d) f(v) for paper (2f|. 



different materials, where the last two data sets clearly 
demonstrate the non-monotonic nature of sliding friction, 
and the first one presumably does not span a sufficiently 
large range of i>'s to detect a minimum. 

Consider now a homogeneous driving stress T d . For 
T d < r m there exists only one stable homogeneous solu- 
tion, the clastic (static) one. Upon increasing r d above 
r m , three solutions exist: the elastic one with v = and 
two steady sliding solutions, one with v < v m (typically 
unstable) and one with v > v m (typically stable). The 
critical point r d = r m corresponds to a bifurcation, which 
suggests a qualitative change in the behavior of the sys- 
tem. At this point we expect steady state propagating 
rupture, in which a solution with v > v m invades an elas- 
tic (static) solution with v — 0, to emerge. Denote the 
propagation velocity of such fronts by c and the one cor- 
responding tO T d =T m by C m in- 

In order to find propagating rupture solutions, and in 
particular to calculate c m i„, we need to couple the fric- 
tion law in Eqs. (JTJ)- (|4]) to an elastic body. It would 
be very difficult to analytically calculate c m i n when the 
body is a 2D medium. Therefore, to gain analytic in- 
sight into the properties of the steady state fronts, we 
assume that the height H of the elastic body (say in the 
y-direction) is much smaller than the spatial scale of vari- 
ation £ of fields along the interface (in the x-direction), 
i.e. we consider a quasi- ID limit. Under these condi- 



tions, integrating the bulk momentum balance equation, 
pdttUi — djUij (where p is the bulk mass density, ct^ 
are the stress components and Ui are the displacement 
components) along the y-direction yields 



Hpd u u(x,t) ~ Hp,d xx u(x,t) 



r(x,t) 



(8) 



where p is the bulk shear modulus and u — u x is the 
intcrfacial shear displacement (slip) that satisfies dtu = v. 
Note that we have omitted constants of order unity in Eq. 
(|8]) and that in the quasi-lD limit both the driving stress 
T d and the friction stress r do not appear as boundary 
conditions, but rather as terms in the "bulk" equation. 

We now look for steady state propagating solutions of 
Eqs. flU and © in which all of the fields take the 
form <?(£), with £ = x — ct, such that a sliding solution 
at £ — ¥ — oo propagates into an elastic solution at £ — ► 
oo. Smoothly connecting these two different solutions 
around £ = provides solvability conditions that allow the 
calculation of c. Since the front propagation transforms 
an elastic region with high contact area to a sliding region 
with lower contact area, these are actually rupture fronts. 

Cmin is being estimated using a scaling calculation in 
which the loading r d is homogeneous and equals to its 
threshold value r m . A self-consistency constraint on the 
quasi-lD formulation is H <C £■, where I is the spatial 
scale characterizing all of the fields in the front solution 
(as denned above). We first use dt = —cd^ to transform 
Eqs. ©, ((3]) and © into the following set of coupled 
ordinary differential equations 

H(£- cpj d t v = r d - t , (9) 
_ cd ^ = a) v _ \A T ei 6 (1 - Tc ) , (10) 



-cd^4> = 1 



M 

D 



(ii) 



We stress that the front velocity c in these equations is 
not a-priori known, but is rather a "nonlinear eigenvalue" 
of this problem. It is determined from the condition that 
the spatially-varying propagating solution properly con- 
verges to the homogeneous sliding solution at £ — > — oo 
and to the homogeneous elastic solution at £— >oo. 

To proceed, we estimate the right-hand-side of 



Eq. dTUJ) as 



% T el I D and the left-hand-side 



as 



/t(c m in), yielding 



(12) 



We neglect inertia (i.e. assume c m i n <C c 8 = 
\J pj p) and estimate the left-hand-side of Eq. (|9]) as 
—H pv-mj l{cmin)c m in- For the right-hand-side, we esti- 
mate t = r el + r ms by its peak value. The maximum 
r el is obtained by solving Eq. (fTU|) with d^r el = 
to get r^ ax - p a DA {l + iPo)/h. 



of 

We estimate that this 
happens roughly half-way through the front, i.e. £ ~ 
—(-{c m i n )/2 and v^v m /2, so the viscous stress is r ms ~ 
rjv* \og(v m /2v*), where we used v m /v* 3> 1. Therefore, 
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we estimate t~ li DA (1 + il;o)/h + tjv* log (i> m /2u*) and 
use r d — r m ~ A [noD/h + rjv* log (w m /v*)] to obtain a 
scaling version of Eq. (j9|). The result, together with Eq. 
(|12p. leads to the following prediction 



(j, H an 



<jD [n D^ /h-riv*log{2)} 



(13) 



where we remind the reader that "00 quantifies the typical 
age of the interface in the elastic regime prior to slip (see 
above). 

Several features of this central result are noteworthy. 
First, c m in is finite and proportional to v m . Second, it 
is independent of inertia, i.e. it does not scale with the 
elastic wave speed c s = y/jjjp [H| . Finally, c m in depends 
on: (i) the properties of the friction law, e.g. on constitu- 
tive parameters such as the viscous-friction coefficient rj 
and the interfacial elastic modulus liq and on microscopic 
geometric quantities such as D and h, (ii) the bulk geom- 
etry through H, (iii) the normal stress as er -1 / 2 and (iv) 
the bulk shear modulus ii. We expect these features to 
remain qualitatively valid independently of the explicit 
form of the friction law and of dimensionality as long as 
steady sliding friction exhibits a non-monotonic behavior 
(cf. Fig. QJ,), as suggested in (l6j . 



IV. AN EXPLICIT EXAMPLE 
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FIG. 2: From top to bottom, «(£), r(£), r ei (£) and for 
a steady state rupture mode propagating from left to right at 
1.2 m/sec <C c s . The parameters for PMMA were used to- 
gether with r d = r m . The spatial scale of variation of the 
fields is denoted by I. The homogeneous states that the 
rupture mode smoothly connects are marked by horizontal 
dashed lines. 



In order to test the prediction in Eq. (|13[) we should 
determine the friction parameters for a specific material. 
We focus on a glassy polymer (PMMA) which was ex- 
tensively characterized in laboratory measurements, in- 
cluding the interfacial elastic response [13] ■ It was used 
in the most conclusive experiments which demonstrated 
the existence of slow rupture modes [l], [f| • 
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TABLE I: Elastic and frictional parameters for PMMA. 

We extracted the material parameters of PMMA from 
various sources. It seems that frictional parameters are 
sensitive to environmental and experimental conditions. 
Somewhat different results for (supposedly) the same ma- 
terial, sometimes by the same research group, were re- 
ported. However, the trends are robust, as well as the 
order of magnitude of the parameters. With this in mind, 
we compiled a list of parameters with which our model 
satisfactorily describes the data from all of the reported 
sources. 

The interfacial elastic response data of Fig. 2 in [l?} 
indicates that han / ' liq is in the fim scale. The slope of the 



ageing data in Fig. 9 of [8| implies P — LioDb/anh — 0.02. 
The "direct effect" measurement in Fig. 4 of [8] implies 
a = r)v*/aH ~ 0.005. v* was estimated as the lowest slip 
rate for which the logarithmic "direct effect" is observed. 
The data presented in Fig. [Tp above are consistent with 
/3— a ~ 0.015. D was determined to be 0.9/im in [l3| and 
0.4/im in [20J. With these constraints, and using known 
values of independently measured parameters such as /i, 
an, p and t c (which is estimated as the yield stress), we 
fitted all of the aforementioned experimental data. The 
parameters are summarized in Table |U In addition, we 
set cr=l MPa, which is similar to the value used in [llllij]. 
H = lO^tm and tpo = 0.6. 

In Fig. [2] we show a steady state rupture solution ob- 
tained by numerically integrating our model equations in 
©-(HU) using these parameters and the threshold driv- 
ing stress r d — r m . The front propagation velocity at 
threshold, c m %n — 1.2m/sec, is more than three orders of 
magnitude smaller than c s « 1, 600m/sec, qualifying it as 
"slow rupture", and £(c m i n ) is on a mm scale, satisfying 
H<^£ as required by self-consistency. 

We test the prediction in Eq. (|13l) by each time vary- 
ing one parameter on the right-hand-side over a wide 
range and comparing the prediction to the numerically 
calculated c m i n . The results are presented in Fig. [3] and 
exhibit excellent agreement between the analytic predic- 
tion and the numerically calculated values of c min . This 
result clearly and directly demonstrates the existence of 
friction-controlled slow rupture in our model. 
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FIG. 3: The numerically calculated c m i„ vs. the analytic 
prediction appearing on the right-hand-side of equation (|13fl . 
which we denote here as c pr . Different symbols correspond to 
variations of different parameters of the model, as denoted in 
the legend. The dashed line serves as a guide to the eye. 



converted into £(c) using the results of Fig. [4}:). The nu- 
merical results clearly confirm the theoretical prediction, 
demonstrating that indeed slow rupture is significantly 
less localized than rupture at elastodynamic velocities. 
Furthermore, the exponential dependence predicted in 
Eq. (fT4")) is quantitatively verified and the slope agrees 
with — l/aa to within a few percent. These predictions 
should be tested experimentally. 




V. THE SPECTRUM OF RUPTURE FRONTS 



The finite velocity scale c m i n implies there are no so- 
lutions with c < c m in ) i-e. the existence of a "forbidden" 
range of velocities in the spectrum of steady state rup- 
ture modes 16]. In Fig. |4ji we show the full spectrum of 
rupture propagation velocities as a function of T d > r m . 
Indeed, there are no solutions with c < c m i n and there 
exists a continuum of states between c m in and the elastic 
wave speed c s . This continuous spectrum qualitatively 
agrees with recent high-resolution laboratory measure- 
ments [l|, reproduced here in Fig. |4p. These measure- 
ments, though not obtained under globally homogeneous 
loading and were done in 2D, directly demonstrate the 
existence of a threshold driving stress, a minimal slow 
rupture velocity and saturation at an elastic wave speed. 
A detailed quantitative comparison to the experiments 
requires fully 2D calculations which are currently under- 
way and will be reported elsewhere. 

Upon increasing r d sufficiently above r m , rupture trav- 
els at a non-negligible fraction of the sound speed and we 
can no longer neglect the inertial term in Eq. The 
sliding velocity in the sliding region of a propagating rup- 
ture solution corresponds to the larger v of the two solu- 
tions of r d = f(v)a, where f(v) is given in Eq. ([5]). For 
H>« m it reads v~v* exp(T d /aa — f3/ab). The analog of 
relation (IT21 for the inertial regime, c~c s ^>c m i n , reads 
£(c r ^c s )^Dc/v, which immediately implies 



,-T d / a < 



(14) 



The strong inequality results from the exponential decay 
of 1(c) with T d in the inertial regime and the typically 
small value of a (~ 0.01, see above). This result predicts 
that slow rupture is much less spatially localized as com- 
pared to ordinary fast rupture. To test this prediction 
we plot in Fig. 0J: £ vs. T d /r m (which can be readily 
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FIG. 4: (a) Semi-logarithmic plot of c/c s vs. r d /r m , for a 
fixed a. The same results are shown in linear-linear scale 
in the inset. The horizontal dashed lines denote c m in and 
c s S> Cmin , and the vertical dashed line denotes the minimal 
driving stress r d — r m . (b) Semi-logarithmic plot of c/cl vs. 
t/ct in the PMMA experiments of [l|] (data courtesy of O. 
Ben-David and J. Fineberg). The same data are shown in 
linear-linear scale in the inset, cl is the longitudinal wave- 
speed and r/cr is rescaled such that the minimal value below 
which no rupture modes were observed equals unity. c m i„ 
in the experiment was of the order of 10 m/sec. (c) Semi- 
logarithmic plot of I vs. r d /r m for the spectrum shown in 
panel (a). 



VI. SUMMARY AND CONCLUSIONS 

Our results, based on a rate-and-state type friction law, 
show that slow rupture is a well-defined and generic state 
of frictional interfaces. The non-monotonic dependence 
of the steady state sliding friction on the slip velocity 
is shown to be directly related to the emergence of a 
new, friction-controlled, velocity scale c m i„ below which 
no steady state rupture can propagate. Furthermore, our 
analysis demonstrates that rupture states span a con- 
tinuum, from friction-controlled slow rupture to inertia- 
limited, earthquake-like, fast rupture [4j. We suspect 
that transient rupture modes observed under complex, 
spatially inhomogcncous, conditions arc short-lived ex- 
citations of these steady rupture states, as was argued 
and demonstrated in a specific context in [l6| . If true, 
steady state rupture fronts play a role analogous to "nor- 
mal modes" or "eigenstates" in other dynamical contexts. 
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The results presented were shown to be in qualitative 
agreement with the recent high-resolution real-time labo- 
ratory measurements of [H[H|. A quantitative comparison 
to the experimental data requires 2D calculations which 
are currently underway. We believe that the developed 
ideas are applicable to geophysical systems, though such 
questions go well beyond the scope of the present work. 
The possible geophysical relevance of our results (e.g. to 
slow/silent earthquakes) should be carefully examined in 
a future investigation. 
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